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Abstract 


The overabundance of the red and massive candidate galaxies observed by the James Webb Space Telescope 
(JWST) implies efficient structure formation or large star formation efficiency at high redshift z~ 10. In the 
scenario of a low or moderate star formation efficiency, because massive neutrinos tend to suppress the growth of 
structure of the universe, the JWST observation tightens the upper bound of the neutrino masses. Assuming A cold 
dark matter cosmology and a star formation efficiency €[0.05, 0.3] (flat prior), we perform joint analyses of Planck 
+JWST and Planck+BAO+JWST, and obtain improved constraints Xm, < 0.196 eV and » 7m, « 0.111 eV at 
95% confidence level, respectively. Based on the above assumptions, the inverted mass ordering, which implies 
Xm, > O.leV, is excluded by Planck+BAO+JWST at 92.7% confidence level. 
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1. Introduction 


The standard hot big bang cosmology predicts a cosmic 
neutrino background (CvB), which decoupled from the thermal 
bath in the early universe at a temperature ~MeV. The 
subsequently redshifted momenta of CvB follow an ultra- 
relativistic Fermi-Dirac distribution with negligible chemical 
potential and an effective temperature 1.95(1 + z)K, where z 
denotes the cosmological redshift. Although direct detection of 
CvB is yet unrealistic, the existence of CvB has been indirectly 
confirmed by the observations of primordial abundances of 
light elements and the cosmic microwave background (CMB). 

The difference between the squares of the three mass 
eigenstates has been measured in neutrino flavor oscillation 
experiments (Fukuda et al. 1998; Abe et al. 2014). These 
experiments constrain the sum of masses of three neutrino 
species, »7m,, to be >0.06eV for the normal mass ordering 
m, ~ m, & mz, and >0.10 eV for the inverted mass ordering 
ma Km, % m». Distinguishing between the two mass-ordering 
scenarios, and hence accurate measurements of the sum of 
neutrino masses, are important for understanding the origin of 
neutrino masses. 

At z zz 1100 where primary CMB anisotropies are generated, 
the background temperature is 72:3000K z0.3eV. If 
Xm, 20.3 eV, neutrino masses will have an impact on the 
primary CMB via the early integrated Sachs Wolfe effect. In 
the late universe, massive neutrinos act as a hot dark matter 
component, which tends to suppress the growth of large scale 
structure at small scales ($100 Mpc) and alter the CMB lensing 
effect. The third generation CMB experiment Planck mission 


constrained the total neutrino masses to be 5» 7m, < 0.241 eV 
(TTTEEE + lowE + lensing) at 95% confidence level (CL) in 
the standard A cold dark matter (ACDM) model (Aghanim 
et al. 2020). Further including baryon acoustic oscillation 
(BAO) data, which breaks the degeneracy between neutrino 
masses and the background expansion history of the universe, 
pushes the upper limit of neutrino masses to ?7n, < 0.121 eV 
(95% CL). This upper bound already puts inverted mass 
ordering under pressure (excluded at 90.2% CL). 

The standard cosmological structure and galaxy formation 
are recently challenged by red and massive candidate galaxies 
from James Webb Space Telescope (JWST). Labbé et al. 
(2023) found six candidate massive galaxies (stellar mass Loe" 
solar masses) at 7.4 < z < 9.1. This finding suggests that either 
the star formation efficiency (SFE) significantly exceeds its 
typical values in low-redshift galaxies, or the halo mass 
function is about 2c higher than the prediction of the standard 
ACDM cosmology (Haslbauer et al. 2022; Wang et al. 2023a; 
Boylan-Kolchin 2023; Lovell et al. 2023; Qin et al. 2023). 
Although at the current stage a large SFE cannot be 
theoretically excluded (Zhang et al. 2022; Qin et al. 2023), a 
large SFE at z> 7 tends to enhance the UV radiation from 
galaxies, causing a faster reionization process that is in tension 
with observation (Goto et al. 2021; Wold et al. 2022; Lin et al. 
2024). For this reason, we will take the more widely accepted 
conservative assumption of SFE 0.1, and discuss the 
cosmological implications. Because massive neutrinos tend to 
suppress the number of massive halos in which massive 
galaxies are born, a large ? 7m, increases the tension between 
the ACDM model and JWST data. Thus, adding JWST data to 
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the cosmological analysis, as the present work aims to do, helps 
to tighten the upper bound of neutrino masses. 

This paper is organized as follows. Section 2 constructs the 
likelihood of massive galaxy counting and describes the 
numerical tool for the likelihood evaluation. In Section 3 we 
apply the likelihood to the Labbé et al. (2023) data and obtain 
Planck--JWST and Planck+BAO+JWST joint constraints on 
>d=m,. Section 4 discusses and concludes. 


2. Likelihood of Massive Galaxy Counting 


The aim of this section is to extend the likelihood of massive 
galaxy counting in Wang et al. (2023a), which we sketch 
below, to massive neutrino models. 

To assess the impact of massive neutrinos on galaxy 
formation, we adopt the extended Press-Schechter ellipsoidal 
collapse model (Press & Schechter 1974; Sheth et al. 2001; 
Sheth & Tormen 2002) to compute the expected abundance of 
dark matter halos. It predicts the halo mass function, which is 
the comoving halo number density per mass interval, to take 
the form 


dN p dinv 
= 1 
dM ima? (1) 


Here, p is the average density of the background and 
v= 6,/o(M, z), where 6,— 1.686 corresponds to the critical 
linear overdensity and o(M, z) is the mass fluctuation at scale 
R = (38M/ 470). The simulation-calibrated f(v) factor is given 
by 


fa) = 2a( 4 + Jes. (2) 


where v’ = Jav, a — 0.707, A — 0.322 and g=0.3 (Sheth & 
Tormen 1999). For massive halos we assume that the fraction 
of baryonic mass in a massive dark matter halo is fp, = Q5/€,,, 
where Q, and Q, are the parameters of baryon and matter 
abundance. The stellar mass M, is then connected to the halo 
mass Mhaio via M. = efbMhaio, Where e is the SFE. This relation 
allows to convert a stellar-mass threshold M. cu; to a halo-mass 
threshold Mhalo,cut = ae Assuming that each massive halo 
has a massive central galaxy, we can write the expected number 
of massive galaxies above the stellar-mass threshold and in a 
selected comoving volume as 


M f we AN WP oe (3) 


Nn) = 47 d 
\ m) They Mhato,cut *min dM dzdQ 


where the selected comoving volume is defined by the redshift 
interval [Zmin, Zmax] and the sky fraction fa. Following Wang 
et al. (2023a) we take Zmin = 7 and Zmax = 10. The survey area 
is 38 arcmin’, which gives fa = 2.56 x 107. ahs comoving 


volume per redshift interval per solid angle is specified 


? de pe 
by the cosmology. For the spatially flat ACDM model, to a 
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Figure 1. From top to bottom are the halo mass functions at z= 8.5 for 
dm, = 0.06 eV, 0.12 eV, and 0.24 eV, respectively. The other cosmological 
parameters (Qm, Qp, Ho, As, Ns) are fixed at Planck best-fit values. 


very good approximation, we have 
dV c? 
dO HSW + 2? + 1 = Om 


2 
dz! 


i EU st zy + 1 m Qin 


(4) 


where c, Qm and Ho are the speed of light, the total matter 
(CDM + baryon + neutrinos) abundance parameter and the 
Hubble constant, respectively. 

At z<10, massive neutrinos become non-relativistic and 
contribute a small (~10~*) fraction to the total energy budget 
of the universe. However, the correction to Equation (4) due to 
neutrino masses, namely the volume effect, is still negligible 
for our analysis. The major impact of neutrino masses is on the 
matter power spectrum, integration of which gives o(M, z) in 
Equation (1). Figure 1 shows the dependence of halo mass 
function, evaluated at z — 8.5 (the center of the redshift bin 
used in the likelihood), on the sum of neutrino masses. 

The stellar mass and the redshift uncertainties propagate to 
the uncertainty in Nobs, the number of observed galaxies above 
the stellar mass threshold and within the selected volume. 
Meanwhile, due to cosmic variance and Poisson shot noise, the 
theoretical prediction of Np, (denoted as MN) is also 
probabilistic. The likelihood of the theory is then characterized 
by the probability of finding No, < Min, which is given by 


P (Nos < Nn) = 5 P (Nos) »» P(Na). (5) 
Novs=0 Na — Nobs 


The distribution function P(N) is based on the cosmological 
model and the star formation efficiency. For the purpose of 
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Figure 2. Left panel: dependence of (Nn) on 27m, and SFE; other cosmological parameters (Qn, Qg, h, As, ns) are fixed at the Planck best-fit values. Right panel: 
P(Nops) from the 13 red and massive candidate galaxies in Labbé et al. (2023) and P (Nn) from the Planck best-fit ACDM cosmology. 


using rare-object statistics, we are only interested in the 
(Nin) S; O(1) case, where Ny, approximately follows a Poisson 
distribution. The expectation-value parameter of the Poisson 
distribution, denoted as A, has an uncertainty due to cosmic 
variance (Trenti & Stiavelli 2008). Marginalizing over the 
cosmic variance of A that is denoted as a, we obtain the 
distribution function of Np, 


A- (Qu)? ANa 


oo 1 — 
P (Na) = f e ñ eA 
(Min) EE Nal 


(6) 


The distribution function P(N ;) is obtained by randomly 
sampling stellar masses and redshifts of the candidate galaxies 
with the summary statistics given in Labbé et al. (2023), and 
counting the number of galaxies in the selected volume. The 
precise values of P(N»s) rely on the functional form of the 
systematic errors of the stellar masses and redshifts of the 
galaxy candidates. We adopt the triangular distribution 
function, which is shown to be a more conservative estimation 
than smooth distributions, e.g., the skew-normal distribution. 
See Wang et al. (2023a) for more detailed description of the 
algorithm and its theoretical explanation. 

We integrate the cosmological linear perturbation Boltzmann 
solver CAMB (Lewis et al. 2000; Howlett et al. 2012) into the 
P(Nops < Mn) calculator in Wang et al. (2023a). As shown in 
the left panel of Figure 2, the expected number of massive 
galaxies (Nm) decreases with the sum of neutrino masses, 
while the other cosmological parameters (Q,, = 0.3158, 
Q, = 0.04939, h=0.6732, A,=2.101 x 107°, n, — 0.9661) 
are fixed at the Planck best-fit values (Aghanim et al. 2020). 
Here the reduced Hubble constant h is defined by 
Ho— 100^ kms! Mpc '; Omn and Q, are the total matter 
(CDM + baryon + neutrinos) and baryon abundances; A, and 
ns are the amplitude and spectral index of the power spectrum 


of the primordial comoving curvature perturbations. In the right 
panel of Figure 2 we show the derived data, namely P(Nops) for 
the 13 candidate galaxies from Labbé et al. (2023), and the 
theoretical P (Np) for two representative values SFE = 0.1 and 
SFE = 0.3, respectively. For the case of SFE — 0.1, there is a 
tension between the theory and the data, as the theory prefers a 
smaller number of massive galaxies (N « 2) than what the data 
indicates (N = 2). For SFE = 0.3, the theory is well consistent 
with the data as P(Nq) is comparable or greater than P (Nobs) 
for all N's. 


3. Results 


For the given set of observed galaxy candidates and selected 
volume (7 <z < 10 and survey area 38 arcmin2), the cosmol- 
ogy-dependent P(N5&, < Mn) in Equation (5) gives a likelihood 
of cosmological parameters and SFE. For more detailed 
information of the data set, such as the redshifts and stellar 
masses of the candidate galaxies, we refer the readers to Wang 
et al. (2023a). 

Due to the degeneracy between cosmological parameters, 
JWST data alone cannot give a meaningful constraint on ? 7n,. 
We therefore analyze Planck+JWST and Planck+BAO 
+JWST jointly by performing importance sampling with the 
Planck and Planck--BAO Monte Carlo Markov chains from the 
Planck Legacy Archive (https:/ /pla.esac.esa.int/pla/). The 
base parameters, on which flat priors are assumed, are the 
baryon density Qr, the CDM density QU), the reionization 
optical depth 7, the angular extension of sound horizon on the 
last scattering surface 0, the sum of neutrino masses ?7m,, the 
primordial power spectrum parameters In A, and n,, and the 
SFE. Before doing cosmological analysis, we have to assume a 
prior on SFE. Unless otherwise stated, we focus our study on a 
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Table 1 
Marginalized 16 Constraints on Cosmological Parameters; for ? 7m, the 95% CL Upper Bounds are Shown 
Planck Planck+JWST Planck+BAO Planck+BAO+JWST 
Xm, [eV] «0241 «0.196 «0.121 «0.111 
Og 0.807 0.016 0.811 + 0.014 0.813 + 0.010 0.815 + 0.009 
Ho [km s~! Mpc !] 67.1 + 1.05 67.3 + 0.92 67.8 + 0.51 67.9 + 0.50 
1000,42 2.237 + 0.015 2.238 + 0.015 2.242 + 0.014 2.242 + 0.014 
Qh 0.1201 + 0.0013 0.1200 + 0.0012 0.1193 + 0.0009 0.1194 + 0.0009 
1008 1.040 82 + 0.00032 1.040 90 + 0.00031 1.041 00 + 0.00031 1.040 94 + 0.00031 
T 0.0547 + 0.0076 0.0553 + 0.0076 0.0553 + 0.0073 0.0559 + 0.0073 
In (1019A,) 3.046 + 0.015 3.048 + 3.045 + 0.015 3.046 + 0.015 
ns 0.9648 + 0.0042 0.9652 + 0.0042 0.9666 + 0.0037 0.9667 + 0.0037 
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Figure 3. Marginalized 68% CL and 95% CL constraints on (?7,, os) (left panel) and (?7n,, Ho) (right panel). The units are eV for ? 7m, and km s! Mpc! for Ho. 


The prior on star formation efficiency is SFE € [0.05, 0.3]. 


low or moderate SFE, namely a flat prior SFE c [0.05, 0.3], 
though other options will also be briefly discussed. 

The constraints on some key parameters (base or derived), 
before and after incorporating JWST likelihood, are shown in 
Table 1 for a comparison. The 95% CL upper bound of ? 7m, is 
improved by 19% (without BAO) and 8% (with BAO), 
respectively. Although the improvement is insignificant, the 
inclusion of JWST data pushes the upper bound to 
Xm, < 0.111 eV, very close to the Xm, > 0.1 eV limit implied 
by the inverted mass ordering. Integrating over the posterior 
distribution of ? 7m,, we find that the inverted mass ordering is 
excluded by Planck+BAO+JWST at 92.7% CL. 

The addition of JWST data also has an impact on os, the root 
mean square of density fluctuation in a top-hat spherical 
window with radius 84^! Mpc. This is because the suppression 
effect due to neutrino masses is degenerate with the amplitude 


of matter density fluctuations. This degeneracy can be clearly 
seen in the left panel of Figure 3, which shows the marginalized 
lo and 2c contours of ?7m, and og. In the right panel of 
Figure 3 we show the well known geometric degeneracy 
between ? 7m, and Ho, which allows the inclusion of BAO (or 
any other geometric probes at low redshift) significantly 
improve the CMB-alone constraint on neutrino masses. The 
figures are made with the publicly available GetDist 
package (Lewis 2019). 

The solid and filled contours in Figure 4 show the 
marginalized 1c, 2c and 30 constraints on SFE and Xm. 
Although seemingly the upper bound of Xm, tightens when 
SFE is lower, it is not the case. This is because the likelihood is 
very sensitive to SFE, and thus the contours vary significantly 
when we change the prior on SFE. For example, when the 
upper bound of SFE is changed to 0.2, the contours become the 
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Figure 4. Marginalized 1c, 2c and 3c constraints on SFE and 57m, with 
Planck+BAO+JWST. The prior on star formation efficiency is SFE c [0.05, 
0.3]. Orange dashed lines show how the contours vary when the upper bound 
of SFE is changed to 0.2. 


Table 2 
Dependence of the 95% CL Upper Bound of 7m, on the Prior on SFE; Data 
Sets: Planck+BAO+JWST; the Last Column Shows the Relative Improvement 
on the 27m, Bound Against the Planck--BAO (Without JWST) Result, Namely 
Em, < 0.121 eV (95% CL) 


9596 CL Upper Bound of Improvement Against 


Prior on SFE Em, Planck+BAO 
SFE € [0.05, 0.1] Xm, < 0.118 eV 2.5% 
SFE € [0.05, 0.2] Xm, < 0.111 eV 8.3% 
SFE € [0.05, 0.3] Xm, < 0.111 eV 8.3% 
SFE € [0.05, 0.4] Ym, < 0.114 eV 5.8% 
SFE € [0.05, 0.5] Xm, < 0.117 eV 3.3% 
SFE = 0.1 dim, < 0.114 eV 5.8% 
SFE = 0.2 dm, < 0.110 eV 9.196 
SFE — 0.3 Mm, « 0.113 eV 6.696 
SFE = 0.4 Mm, < 0.118 eV 2.5% 
SFE = 0.5 Ym, < 0.121 eV 0% 


orange dashed lines in the figure. When we marginalize over 
SFE, the upper bound of ? 7m, does not change much. See also 
Table 2 for concrete numbers. 

By varying the prior on SFE, we also find an interesting 
phenomenon that, as Table 2 shows, the best improvement of 
Xm, bound is achieved around SFE ~ 0.2. Physically this is 
because a very small SFE makes the likelihood approach a 
cosmology-independent constant P — P(Nə»s = 0), that all 
candidate galaxies are interpreted as false detections due to 
systematic errors in redshift and mass, and a very large SFE 
makes the likelihood approach another cosmology-independent 
constant (P — 1), that theoretical prediction significantly 
exceeds the observed number of galaxies. 
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4. Conclusions and Discussion 


In this work, we incorporate massive neutrinos into the 
numerical tool in Wang et al. (20232) to evaluate the likelihood 
of massive galaxy counting. At the current stage, the JWST 
candidate galaxies alone cannot constrain cosmology. How- 
ever, in joint analyses with Planck or Planck+BAO, the 
addition of JWST massive galaxy counting improves the upper 
bound of neutrino masses by a few percent and puts the 
inverted mass ordering under more pressure. 

While the Planck+BAO+JWST constraint gives the most 
stringent upper bound on neutrino masses, the results are yet 
very model-dependent. In particular, our results rely on SFE 
prior at z ~ 7 , namely SFE < 0.3. If a larger SFE > 0.4, which 
is unlikely in the usual scenario of galaxy formation and 
evolution, is assumed, the likelihood becomes insensitive to 
cosmology and can no longer improve the 57m, constraint. 
However, there are alternative scenarios where an enhanced 
SFE at high redshift does not accelerate the reionization 
process, and therefore it is consistent with observations. For 
example, Lin et al. (2024) considered large SFE (20.4) at high 
redshift and warm dark matter in the keV mass range. In this 
picture, the warm dark matter has no impact on the abundance 
of very massive halos, and only suppresses the low-mass end of 
halo mass function. In other words, the keV warm dark matter 
reduces the number of low-mass galaxies that emit UV 
photons. This effect cancels the enhancement of UV radiation 
from each galaxy due to a large SFE, and keeps the cosmic 
reionization history in line with observations. There are also 
other options, such as fuzzy dark matter, which can suppress 
the abundance of small halos and galaxies and lead to the same 
phenomenon (Gong et al. 2023). 

Other cosmological probes, the type Ia supernovae, local 
distance-ladder measurements of Hp, weak gravitational 
lensing, redshift-space distortion, etc., may also tighten the 
current upper bound of neutrino masses. See e.g., Aghanim 
et al. (2020) and discussion therein. Ideally, one would like to 
join all the available data and pushes the bound ? 7n, to the best 
limit. However, at the current stage, these cosmological probes, 
including the JWST high-z massive candidate galaxies, are not 
all in concordance and may suffer from unknown 
systematics (Asgari et al. 2021; Riess et al. 2021; Huang 
et al. 2022; Wang et al. 2023b). It is thus not so favorable to 
add too many cosmological data into the pool at once, and 
caution should be taken for the cosmological bounds on 
neutrino masses. Nevertheless, the growing JWST data bring a 
hope to cross check the constraints on neutrino masses in very 
different perspectives (tail statistics, perturbation-based, non- 
linear physics). 
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